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COMPUTER  PROGRAM  FOR  CALCULATING  DIELECTRIC  PROPERTIES 
OF  LOW-  OR  HIGH- LOSS  MATERIALS  FROM  SHORT-CIRCUITED 
WAVEGUIDE  MEASUREMENTS1 


S.  0.  Nelson,2  C.  W.  Schlaphoff,3  and  L.  E.  Stetson2 


A  general  computer  program  with  a  number  of  unique  features  was  developed 
for  precise  calculation  of  dielectric  properties  of  materials  from  measure- 
ments obtained  on  short-circuited  coaxial,  cylindrical,  or  rectangular  wave- 
guides.    Dielectric  properties  of  either  high-  or  low-loss  materials  are  cor- 
rectly determined  from  standing-wave  ratios  (measured  with  and  without  the 
sample  in  the  shorted  end  of  the  waveguide)  and  the  shift  of  the  voltage 
standing-wave  minimum  (or  the  current  standing-wave  maximum)  resulting  from 
the  presence  of  the  sample.     Corrections  are  included  for  waveguide  wall 
losses  and  for  the  influence  of  the  slot  on  the  guide  propagation  constants 
for  rectangular  and  cylindrical  waveguide  slotted  sections.     Important  prin- 
ciples of  the  measurement  and  calculation  methods  are  outlined.  Graphical 
illustrations  and  necessary  mathematical  relationships  are  presented.  A 
listing  of  the  program  is  included  along  with  concise  examples  of  input  and 
output  data. 

INTRODUCTION 

The  dielectric  properties  of  materials  are  important  for  numerous  ap- 
plications.    In  the  study  for  which  this  computer  program  was  developed, 
data  were  needed  over  a  wide  range  of  frequencies  on  the  dielectric  proper- 
ties of  grain  and  insects  to  identify  the  most  promising  frequency  ranges 
for  possible  insect  control  by  selective  dielectric  heating  of  the  insects 
in  infested  grain  (8).^    The  program  and  methods  of  measurement  described 
here,  however,  are  applicable  to  a  wide  range  of  materials.     They  are  also 
useful  over  a  wide  range  of  frequencies  where  slotted-coaxial-line  and 
slotted-waveguide  measurement  techniques  are  applicable,  approximately  100 
MHz  to  30  GHz  or  higher. 

The  dielectric  properties  of  interest  here  are  the  permittivity  or 
dielectric  constant,  the  dielectric  loss  factor,  the  loss  tangent  of  the 
dielectric,  and  the  a-c  conductivity.     The  reader  is  referred  to  other  sources 
for  detailed  conceptual  information  and  definitions  based  on  circuit  concepts 
(4,  6)  or  electromagnetic  field  concepts  (3,  7,  8).     The  complex  relative 
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permittivity  (complex  relative  dielectric  constant)  is  represented  as 

er  =  Ep  -  j£p  =  -  j  tan  6),  where  £p  is  the  real  part  of  the  complex 

permittivity,  frequently  called  the  dielectric  constant;  the  imaginary  part, 

,   is  the  dielectric  loss  factor;  and  tan  6  =  e^'/ is  known  as  the  loss 
tangent  (also  called  the  dissipation  factor).     The  a-c  conductivity, 
a  =  oje  e^',  and  e^J  include  all  energy-dissipating  mechanisms  observed  in  the 
dielectric  as  a  result  of  its  subjection  to  the  forces  of  the  alternating 
electric  field.     Here  eq  is  the  permittivity  of  free  space,  8.85^+  x  10-^2 
farad/m.  ,  and  to  represents  the  angular  frequency,  2-nf,  where  f  is  frequency 
in  Hz . 

By  using  a  measurement  system  with  a  slotted-line  or  slotted-wave guide 
section,  it  is  possible  to  obtain  information  from  which  the  dielectric 
properties  of  materials  may  be  calculated.     Roberts  and  von  Hippel  (10) 
originally  reported  the  method  used.     The  voltage  standing-wave  ratio  (VSWR) 
in  front  of  a  sample  material,  placed  at  the  end  of  a  coaxial  line  or  wave- 
guide terminated  in  a  short  circuit,  is  determined  with  the  sample  in  place. 
The  VSWR  is  also  measured  with  the  sample  removed.     The  shift  in  the  posi- 
tion of  the  voltage  minimum  (node)  due  to  the  presence  of  the  sample  is  also 
determined.     This  information,  together  with  the  sample  dimensions,  frequen- 
cy, and  characteristics  of  the  waveguide,  permits  determination  of  the  di- 
electric properties.     Solution  of  a  complex  transcendental  equation  required 
the  use  of  charts  of  complex  hyperbolic  functions  or  a  series  approximation. 
Dakin  and  Works  (2)  reported  simplified  calculation  procedures  for  this 
method  based  on  approximations  valid  for  materials  when  tan  6  <  0.1.  They 
also  included  a  correction  for  the  wall  losses  in  the  waveguide  in  deter- 
mining the  loss  tangent  of  the  sample  material.     Westphal  (12),  who  included 
a  more  precise  procedure  for  correcting  for  wall  losses  in  the  waveguide, 
described  further  refinements  for  the  method.     A  complete  description  of 
the  principles  of  this  method  and  methods  of  measurement  and  calculation, 
along  with  derivations  of  important  equations ,  has  been  assembled  in  one 
reference  (8).     Only  those  relationships  necessary  for  understanding  and 
describing  the  computer  program  will  be  presented  here. 

PRINCIPLES  OF  THE  MEASUREMENT  METHOD 

Fig-.  1  represents  the  short-circuited  end  of  a  rectangular  waveguide 
containing  a  sample  of  length  d  in  contact  with  the  short  circuit.     The  in- 
cident wave  traveling  in  the  positive  z  direction  and  the  reflected  wave 
combine  to  form  standing  waves  in  the  waveguide.     A  voltage  standing-wave 
pattern  is  illustrated  in  fig.  2,  where  it  may  be  noted  that  the  wavelength 
in  the  sample,  medium  2,  is  shorter  than  that  in  the  empty  part  of  the  wave- 
guide, medium  1.     Voltage  minima  or  nodes  are  spaced  at  half-wavelength  in- 
tervals in  the  waveguide.     A  more  general  voltage  standing-wave  pattern  is 
illustrated  in  fig.   3,  where  a  comparison  of  the  standing-wave  pattern  with 
and  without  the  sample  in  the  end  of  the  waveguide  is  shown.     In  these  fig- 
ures,       represents  the  guide  wavelength  in  the  empty  waveguide,  hereinafter 
referred  to  as  \g  for  conventional  reasons ,  zs  represents  the  position  of 
the  voltage  node  with  the  sample  in  the  waveguide  (the  sample  node),  and  za 
represents  the  node  position  for  the  empty  waveguide  (the  air  node).  hza 
and  Asg  represent  widths  of  the  air  and  sample  nodes,  respectively,  measured 


Figure  2. — Relationship  of  voltage  standing-wave  pattern  to  shorted 
waveguide  with  dielectric  sample. 


_  1+  _ 


Figure  3. — Voltage  standing-wave  pattern  relationships  in  empty 
waveguide  and  in  waveguide  with  dielectric  sample  in  place. 

at  specified  levels  above  the  voltage  minimum,  and  na  and  ng  are  integers. 
The  distance  from  the  air-sample  interface  to  the  first  voltage  minimum  out- 
side the  sample  is  designated  as  z  . 

If,  instead  of  a  slotted  line  with  a  capacitively  coupled  probe,  an  in- 
ductively coupled  loop  in  the  plane  of  a  movable  short  circuit  is  used,  the 
current  standing-wave  pattern  may  be  observed.     In  this  case,  the  locations 
of  the  current  maxima  are  measured  and  the  width  of  the  current  peaks  at  a 
specified  level  below  the  maximum  is  used  to  determine  the  standing-wave 
ratio . 
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As  already  mentioned,  the  measurement  of  za  and  zQ  ,  and  the  standing- 
wave  ratios  with  and  without  the  sample  in  place  are  sufficient,  together 
with  known  values  for  d  and  frequency,  or  guide  wavelength,  to  calculate 
the  dielectric  properties  of  the  sample.     Standing-wave  ratios  are  obtained 
from  the  measured  As  values  and  the  dB  levels  above  the  minimum  where  the 
node  widths  are  measured  in  accordance  with  the  following  equation: 
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where  Pr  represents  the  power  ratio  corresponding  to  the  dB  level  employed. 
It  has  been  customary  in  these  types  of  measurements  to  use  3-dB  levels 
(twice-minimum-power  points)  for  simplicity,  since  eq.  1  then  reduces  to 
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The  value  for  A    needed  in  eqs.  Tor  2  may  be  obtained  by  measurement 
of  Xg/2  with  the  slotted  section  or  calculated  from  the  frequency  if  the 
frequency  is  accurately  known, 
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where  A0  is  the  free-space  wavelength  corresponding  to  f  =  <?/A0  (a  =  velocity 
of  propagation  in  free  space,  2.997925  x  108  m./sec),  and  /    and  Ac  are  the 
cutoff  frequency  and  cutoff  wavelength,  respectively,  for  the  waveguide 
(/CAC  =  o). 


For  precise  measurements ,  the  influence  of  the  slot  in  a  slotted  wave- 
guide section  must  also  be  taken  into  account.     Distances  measured  with  a 
slotted  section  will  be  slightly  too  long  because  reflections  from  the  slot 
influence  the  propagation  constant  in  the  guide.     Corrections  for  this  can 
be  made  using  the  following  relationships.     For  rectangular  waveguide  propa- 
gating energy  in  the  TE1Q  mode  (5,  9), 
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where  A^  represents  the  measured  guide  wavelength,  A^  the  true  guide  wave- 
length in  a  nonslotted  section  of  the  waveguide,  w  the  width  of  the  slot,  and 
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a  and  b  the  inside  dimensions  of  the  rectangular  waveguide  (fig.  l).  For 
cylindrical  waveguide  operating  in  the  TE^  mode,5 
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where  T  is  the  inside  radius  of  the  guide  and  t  is  the  first  root  providing 
a  zero  value  for  the  first  derivative  of  the  first-order  Bessel  function  of 
the  first  kind  and  has  the  value  I.8U118. 


The  original  equations  presented  "by  Roherts  and  von  Hippel  (10)  may  be 
written  as  follows : 
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where  e     is  the  complex  dielectric  constant  of  the  sample  and  Y2  =  a2  +  «7^2 
is  the  complex  propagation  constant  in  the  sample,  where  a  and  (3  are  the 
attenuation  and  phase  constants. 


Values  for  Y2^  must  be  extracted  from  eq_.  7  for  use  in  eq.  6.  Deter- 
mination of  the  inverse  standing-wave  ratio,  Xg ,  zQ ,  and  <i,  the  length  of 
the  sample 4  permits  the  right-hand  side  of  eq.  7  to  be  reduced  to  a  complex 
number,  CeJ^,  and  y2^  may  be  obtained  from  charts  of  the  complex  hyperbolic 
function  (10,  12)  or,  in  some  cases,  calculated  from  a  series  approximation 
(10). 


The  simplified  equations,  for  the  case  when  tan  5  <  0.1,  presented  by 
Dakin  and  Works   (2)  may  be  written 
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5Personal  communication,  Howard  E.  Bussey,  National  Bureau  of  Standards, 
U.S.  Department  of  Commerce,  Boulder,  Colorado,  1971. 
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Using  these  approximations,  values  for  $2^  can  ^e  obtained  more  easily  than 
by  reference  to  tables  for  the  function  (tan  X)/X.     The  dielectric  con- 
stant is  then  obtained  using  eq.  8  and  tan  6  is  given  by  another  expression. 
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The  wall-loss  correction  suggested  by  Dakin  and  Works  was  obtained  by  apply- 
ing eq.  10  to  the  empty  waveguide  giving 
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which  is  subtracted  from  the  tan  6  value  obtained  from  eq.  10  for  the  sample 
measurement  to  obtain  the  corrected  tan  6  value  for  the  sample. 

A  more  precise  correction  for  wall  losses  was  described  by  Westphal  (12) 
taking  into  account  differences  in  wall  losses  in  the  sample-filled  part  of 
the  guide  with  the  sample  present  and  with  the  sample  absent.     To  correct  for 
wall  losses  between  the  sample  face  and  the  slotted  section  probe,  the  value 
of  A3S  can  be  corrected  before  using  it  in  eq.  10,  or  in  eq.  1  for  determina- 
tion of  inverse  standing-wave  ratios  for  use  in  eq.  7.     The  appropriate  ex- 
pression, referring  to  fig.  3,  is 


(A3   )  ,    ,  =  (As  )  j  —  ~j —  As 
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12] 


The  wall  losses  in  the  sample-filled  part  of  the  guide  are  then  determined 
from  the  following  relationship  applicable  to  nonmagnetic  materials , 
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where,  for  rectangular  waveguide,  Cg  =  a/ 2b,  and,  for  cylindrical  waveguide, 
Cg  =  0.U2.     This  term  then,  tan  &ws ,  is  subtracted  from  the  tan  6  value  ob- 
tained from  eq.  10  or  from  tan  6  =  e£/ 'e'  using  values  obtained  from  eq.  6. 
In  the  latter  case,  a  corrected  e"  value  is  then  calculated  as        tan  6  using 
the  corrected  value  for  tan  6  of  the  sample.     Losses  in  the  short-circuit 
termination  are  usually  neglected,  but  may  need  to  be  taken  into  account  for 
measurements  on  extremely  low-loss  materials.     Since  these  losses  are  diffi- 
cult to  describe  generally,  no  provision  for  them  is  included  in  the  program. 
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All  of  the  foregoing  equations,  except  eqs.  k  and  5,  are  also  applicable 
to  coaxial  waveguide  where  AQAC  becomes  zero,  since  the  cutoff  wavelength 
is  infinite.     The  slot  in  a  coaxial  slotted  section  has  no  effect  on  the 
measured  guide  wavelength,  so  no  correction  for  its  influence  is  required. 

The  determination  of  zQ  needed  for  use  in  eq.  7  or  in  eqs.  9  and  10  is 
found  from  fig.  3  to  he 

«o  =  K  ~  ns)  V2  '  [d  +  {Za  ~  Zs)]  '  ™ 

Another  equation  from  the  paper  of  Dakin  and  Works   (2)  is  used  for  a 
first  approximation  in  the  computer  program, 
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For  high-loss  materials,  equations  for  the  general  case,  eqs.  6  and  7, 
must  he  used  to  obtain  valid  data.     Since  the  calculations  and  use  of  com- 
plex function  charts  are  tedious,  time-consuming,  and  sometimes  inaccurate, 
a  computer  program  was  developed  to  handle  the  necessary  computations. 

COMPUTER  PROGRAM  DESCRIPTION 

Main  Program 


The  program  was  designed  for  general  use  and  accommodates  input  data 
taken  from  measurements  on  coaxial  lines,  or  cylindrical  or  rectangular  wave- 
guides.    It  properly  performs  the  necessary  calculations  for  either  low-loss 
or  high-loss  dielectric  samples  and  prints  out  the  final  values  for  e^,  e^, 
tan  6 ,  and  a . 

Calculations  were  programmed  for  computation  on  an  IBM  360  Model  65  com- 
puter using _F0RTRAN  IV  programming  language  and  requiring  lUOK  storage  for 
the  GO  step.     (if  necessary,  storage  requirements  could  be  reduced.)  The 
program  listing  is  shown  in  the  appendix.     A  laboratory  data  sheet  with  en- 
tries for  measurements  is  also  included  in  the  appendix.     Input  data  are 
punched  on  three  cards  (fig.  h) .     An  initial  card  indicates  the  number  of 
data  sets  to  be  run  ^(lijnit^oJ^Oj,-- 

A  parameter/ card  for  each  data  set  carries  the  number  of  samples  for  the 
data  s  e tyf4 0 - char  ac  t  e r  sample  description  for  the  data  set,  waveguide  dimen- 
sions, a  reference  dimension,  R,  which  corresponds  to  an  approximate  distance 
between  the  short-circuit  termination  and  the  slotted-line  or  slotted-section 
probe  when  positioned  at  zero  on  its  scale,  the  width  of  the  slot  in  the  slotted 
section  for  use  in  correcting  slotted-section  wavelength  data,  and  instructions 
for  the  detail  of  printout  desired.     For  rectangular  waveguides,  the  horizontal 


*When  compiled  under  IBM  OS/MVT,  the  program  executes  in  hOK  of  core.  Typ- 
ical CPU  time  required  for  processing  50  sample  measurements  is  less  than  20 
seconds.     There  are  no  files  other  than  the  input  (reader)  and  the  output 
(printer)  files. 
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SAMP  L  E        DATA  CARD 


SAMID 
(A) 
0  0  0  0  0  0  0  0 

'  3  3  i  0  6   7  B 


DATE      M . C 


(A) 
0  0  0  0  0  0  0  0 


(A) 
0  0  0  0 


TEMPF  FREQ 

'  (A) 
000000  00  000 


ANODP 
0  0  0  0*0  0  0 


9  7C  71  2?  73  24  75  76  77  78  7*  70  .1  17  77  J*  7',  76  '7  7b  29  (0  47  47  4  7  44  "j  46  47  48  49  70 


ANODM 
0  0  0  0*0  0  0 


SNODP 
0  0  0  0*0  0  0 


SNODM 
0  5  0  0*0  0  0 


DDIMS 
0  0  0  0*0  0  0 


56  57  5a  5?  60  61  67  67  64  65  66  t '  76  69  73  ,1  72  77  74  75 


ADB 
0  0  0  0' 


SDB  EST 
DC 


0  0  0  0 


0  0  0 


\ 


0  0 


PARAMETER  CARD 


JJ 
0  0  G 


SAMPL 
(A) 

0  0  0  0  0  C  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  2  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0 


2  7  4  5  6    1  «  9  I 


79  »?  712  21  77  72  24  27  76  27  28  29  30  .:  27  77  2  1  .5  77  J)  38  29  43  41  42  ' 


tHDIMW 

o  o  a  0*0  o 


VDIMW 

o  o  o  0*0  n 


RDIMW 
U  0  0  0*0  0 


SLOT 
0  0  C  0*0  0 


ALCALC 
0  0  0  0°0  0 


7  f  0  61  6<  67  64  67  66  67  63  77  "0  /'  li  "3  74  72  76  7 


0  0  0  0  0  0  0 


\ 


INITIAL  COD 


LL 


0  0  0  0  3 

:   2  3  4  5 


U  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0  0 

R  C  10  11  17  17  4  *|5  16  17  IE  19  7(1  21  22  23  24  25  26  27  23  29  20  .!  32  33  34  3'-  jC  27  38  it  40  41  42  43  44  45  46  47  48  49  50  51  52  53  54  55  55  57  58  58  60  ci  62  63  64  6'  -76  £7  63  69  70  7'  72  73  74  75  76  77  75  73  80 


Figure  k. — Three  types  of  input  data  cards  for  computer  program,  initial 
card,  parameter  card,  and  sample  data  card,  showing  division  of  fields 
and  automatic  decimal  point  locations. 

and  vertical  dimensions,  a  and  b  (a  >  b) ,  are  entered  on  the  card.  If  the 
diameter  of  a  cylindrical  waveguide  is  entered  in  the  "a-dimension"  column 
only,  the  calculations  are  carried  out  for  the  cylindrical-waveguide  case. 
If  nothing  is  punched  for  the  a  and  b  .dimensions ,  the  calculation  proceeds 
for  the  coaxial-line  case. 


The  third  card  carries  the  individual  sample  data  including  sample 
identification,  date,  moisture  content,  temperature,  frequency,  twice- 
minimum-power-point  readings  [or  readings  taken  at  other  power-ratio  levels 
for  za  and  z  ,  the  air  (empty-waveguide)  and  sample  voltage  nodes],  the  sam- 
ple length,  a,  the  dB  levels  used  for  the  node  width  readings  if  other  than 
3-dB  levels,  and  an  estimate  of  the  value  of  e^,  the  dielectric  constant. 
All  length  measurements  are  in  mm.,  and  f  is  in  GHz. 

Provisions  are  made  in  the  program  for  calculating  frequency  from  ad- 
jacent node  data  taken  on  the  slotted  lines  or  slotted  sections  if  this  is 
more  convenient  or  results  in  a  more  accurate  frequency  determination  for 
any  particular  system.     If  frequency  is  to  "be  determined  in  this  way,  an 
extra  card  is  punched  in  accordance  with  the  sample  data  card  format  and 
placed  in  front  of  the  sample  data  cards  for  which  this  frequency  applies. 
Wo  entry  is  made  in  the  frequency  columns  on  this  card,  "but  data  for  the  two 
adjacent  air  nodes  are  punched  in  the  air-node  and  sample-node  columns.  Wo 
other  entries  are  made  on  this  card.     The  program  then  checks  for  the  pres- 
ence of  data  in  the  frequency  and  sample-length  columns,  and,  if  none  is 
found,  the  calculation  of  frequency  from  the  guide-wavelength  data  is  per- 
formed.    The  check  for  sample-length  data  is  necessary  since  data  cards 
immediately  following  this  card  have  no  frequency  data  punched  either. 

For  coaxial  lines,  frequency  is  simply  determined  as  f  =  299 . 6966/ , 
since  \Q  =  X   ,  and  is  the  distance  "between  adjacent  nodes  in  the 
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air-filled  line.     The  constant  used  here  represents  the  velocity  for  electro- 
magnetic propagation  in  unbounded  air  space.    Using  e^,  for  air  =  1.0006U,6 
the  velocity  of  light  in  air  =  c//e~  =  2.996966  x  108  m./sec,  where  c  is  the 
velocity  of  light  in  vacuum  (2.997925  x  108  m./sec).    Use  of  the  value  for 
velocity  in  air  provides  higher  accuracy  in  the  measurement  results. 

For  either  cylindrical  or  rectangular  waveguide,  the  relationship  from 
eq..  3,  AQ  =  +  (AgAc)2,  must  be  employed  to  find  the  free-space  wave- 

length, AQ.  Before  calculating  AQ,  however,  the  proper  value  for  \g  must 
be  obtained  from  eqs.  U  or  5 ,  depending  upon  the  type  of  waveguide,  which 
correct  for  the  error  in  measured  guide  wavelength  due  to  the  influence  of 
the  slot  in  the  slotted  section.7  Finding  Xg  in  eqs.  It  or  5  requires  the 
solution  to  a  cubic  equation,  so  an  iterative  solution  using  Newton's  ap- 
proximation , 

x    A  ,  =  x    -f(x  )/f'(x  )   ,  [16] 
n  +  1       n     J    n    J      n  ' 

is  employed  in  the  program.     Once  AQ  is  found,  f  is  calculated  as  299.6966/AQ 
in  the  same  way  as  in  the  coaxial  case. 

Calculations  then  proceed,  using  f,  whether  entered  directly  on  the  data 
card,  or  calculated  as  just  explained  using  adjacent  node  measurement  data, 
as  follows : 

Wavelength  in  air,  AQ  =  299.6966//".    Waveguide  cutoff  wavelength,  Ac  =  2a 
for  rectangular  guide  or  Ac  =  1.706293  D,  where  D  represents  diameter  of  cy- 
lindrical guide.     Guide  wavelength,  Xg  =  XQ  Hi  -  (AQ/AC)Z  .     Next,  for  the 
waveguide  case,  a  correction  factor  for  the  influence  of  the  slot  on  2-dimension 
data  taken  on  a  slotted  section  is  calculated  using  eqs.  h  or  5 •     The  air-  and 
sample-node  positions,  z„  and  so5  with  reference  to  the  slotted-line  or 


6The  value  of  e    for  air  at  20°  C.  and  760  mm.  Hg  pressure  with  a  dew- 
point  temperature  of  9°  C.   (l).    The  equation  for  zv  used  by  Brady  (1)  was 
originally  produced  by  Smith  and  Weintraub  (ll) ,  since  ep  =  n2 ,  where  n  is 
the  index  of  refraction.     Their  equation  was  n  =  1^ —  Jp  +  U8l0  |r  j ,  where 
T  is  the  absolute  temperature  in  degrees  K,  p  is  the  total  pressure  in  mil- 
libars, and  e  is  the  partial  pressure  of  water  vapor  in  millibars.     The  ele- 
vation at  Lincoln,  Nebr. ,  is  about  1,200  feet  above  sea  level,  and  the  cor- 
responding mean  atmospheric  pressure  is  730  mm.  Hg.    Under  these  conditions, 
at  UO-percent  relative  humidity  and  76°  F. ,  the  equation  gives  n2  =  er  = 
1.00061.     Since  the  difference  between  /1.0006U  and  /l.0006l  is  less  than 
2  x  10~5,  no  correction  for  the  1,200-ft.  difference  in  elevation  or  8°-F. 
difference  in  temperature  was  justified  in  this  study. 

7If  slot  dimensions  are  not  entered  on  the  parameter  card,  the  program 
automatically  inserts  values  for  the  equipment  used  in  this  study,  3.2  mm. 
for  rectangular  guide  and  1.2k  mm.  for  cylindrical  guide. 
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slotted-section  scale  zero  are  calculated  by  averaging  the  double-power- 
point  readings  (or  readings  taken  at  other  power-ratio  levels),  and  As„  and 
A3S  are  obtained  as  the  absolute  value  of  their  differences ,  so  that  the 
order  of  double-power-point  entries  for  each  node  on  the  data  cards  is  of 
no  importance.     hza  and  hzQ  values  are  then  corrected  for  the  influence  of 
the  slotted-section  slot.     The  node  shift  due  to  the  sample,  za  -  zg,  is 
then  calculated  and  corrected  for  the  influence  of  the  slot.     At  this  point 
in  the  program,  dQ  of  fig.   3,  which  is  needed  later  for  line  and  wall-loss 
corrections,  is  estimated  by  adding  the  reference  dimension  R  (the  approxi- 
mate distance  from  the  short-circuit  termination  to  the  slotted-line  or 
slotted-section  probe  when  set  to  zero  on  its  scale)  to  the  za  value 
referenced  to  the  scale  zero.     Since  dQ  is  an  integral  number  of  half- 
wavelengths  in  the  guide,  the  precise  value  for  d0  is  found  in  the  program 
by  successively  adding  half-wavelengths  until  the  difference  between  n  Xg/2 
and  the  estimate  for  dQ  is  less  than  \g/4.     With  this  determination  of  n, 
d    is  set  equal  to  n  Xn/2. 

Next,  the  provision  for  using  other  than  3-dB  points   (twice -minimum- 
power  points)  for  the  node-width  measurement  is  considered.     If  no  entry  is 
punched  in  the  dB-level  columns  of  the  data  card,  the  program  automatically 
inserts  3  dB.     The  usual  procedure  is  to  use  3-dB  levels  to  simplify  calcu- 
lations.    However,  better  accuracy  may  be  provided  by  using  levels  greater 
than  3  dB  for  very  high  standing-wave  ratios,  and,  in  practice,  some  samples 
were  found  to  provide  such  a  good  impedance  match  that  3-dB  levels  could  not 
be  used  because  E^n  was  less  than  3  dB  below  5^     .     Therefore,  the  calcula- 
tions are  programmed  to  provide  for  any  possible  dB  level  in  measuring  Asa 
and  A20.     Since  a  greater  dB  level  will  always  be  possible  for  the  air-node 
measurement,  the  program  converts  the  Asa  value  to  the  proper  value  for  the 
dB  level  used  for  the  measurement  of  hzs  if  they  are  measured  at  different 
dB  levels.     Eq.  1  provides  a  value  for  the  VSWR  for  any  dB  level  or  power 
ratio.     The  dB  level  is  converted  to  power  ratio,  Pp  =  [E{ z  ) /Em^n] 2 ,  as  fol- 
lows for  use  in  eg.  1:     dB  =  10  log10  Pp  =  {10  In  Pp) / {in  10),  so 
Pr  =  e° ' 2302585 ( dB ) _     After  finding  the  VSWR  value  from  eq.  1,  the  hza  value 
corresponding  to  the  dB  level  at  which  hzQ  was  measured  is  obtained  by  solv- 
ing eq.  1  for  hza  using  Pr  for  the  sample  dB  level. 

Since  Asg  must  be  corrected  for  guide  or  line  losses  between  the  sam- 
ple and  the  probe,  the  ratio  d'/d&  (see  fig.  3)  is  calculated  and  used  in  eq. 
12  to  provide  the  value  of  A3_  corrected  for  these  losses.     The  inverse 

o 

standing-wave  ratio  with  the  sample  in  the  guide  is  then  calculated  using  eq. 
1. 

Next,  the  program  calculates  zQ  according  to  eq.  ih  by  successively  add- 
ing half-wavelengths  to  -[d  +  {za  -  zg ) ] ,  if  it  is  not  already  positive,  until 
the  result  is  positive.  If  a  sample-node  position  is  used  which  is  farther 
away  from  the  short  circuit  than  the  air-node  position,  {za  -  zg)  will  have 
a  negative  value,  and  if  | za  -  zg |  exceeds  d,  the  quantity  ~{d  +  {&a  -  3g)] 
will  be  positive.  This  value  is  then  taken  for  zQ  in  the  program.  Even  if 
a  sample  node  should  be  selected  which  is  more  than  beyond  the  air  node, 

the  value  obtained  in  this  way  for  zQ  will  work  all  right  in  the  program, 
since  tan  (  2tt  zq/X   )  is  equivalent  to  tan  [2u{zQ  +  \  /2)/\q]  as  a  consequence 
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of  the  identity,  tan  (0  +  tt)  =  tan  0.     Thus,  any  air  node  and  any  sample 
node,  without  reference  to  relative  position  or  distance  apart,  will  suf- 
fice for  the  measurement. 

In  order  to  use  eq.  T  for  computation,  the  right-hand  side  must  "be  ra- 
tionalized and  separated  into  its  real  and  imaginary  components.     This  may 
"be  rewritten 

tanh  y2d          .        X    tan  U($2  -  1)  -  J  X    $(l  +  tan2  V)             X  {A  -  jB) 
 =  qqO  £  =  _J2  2  =   2.  

^2d  2itd{l  +  4>2  tan2  U)  2ird{l  +  (j>2  tan2  U) 

[IT] 

where  tf>  represents  the  inverse  standing-wave  ratio  and  U  represents 
(  2tt  sq)/X~.     For  computer  solution,  close  estimates  for  the  components  of 
y2d  =  a2a  +  jS>2d  are  needed  and  then  they  can  be  adjusted  until  eq.  IT  is 
satisfied.     For  given  values  of  a2d  and  32<i,  an  expansion  of  (tanh  y2d)/y2d 
is  used  to  calculate  real  and  imaginary  parts  in  a  subroutine  called  ACOMP. 

An  initial  estimate  of  32<5  is  obtained  from  eq.  8,  which  is  valid  for 
the  low-loss  case,  by  solving  for  32<^  using  the  estimated  value  from  the 
data  input  for  e'.     In  the  program,  if  no  value  is  entered  on  the  parameter 
card  for  the  estimated  dielectric  constant,  a  value  of  2.5  is  automatically 
taken.     Provision  is  also  made  in  the  program  to  use  the  dielectric  constant 
calculated  for  the  previous  sample  as  the  estimated  dielectric  constant  for 
the  next  sample  in  a  series  of  similar  samples  if  no  estimated  value  is 
punched  on  the  data  .card.     Since,  in  the  lossless  case,  ct2  =  0  =  £ , 
(tanh  y2d)/y2d  =  Ce^^  becomes,  in  accord  with  eq.  9,  (tan  82d)/$2d  =  C.  C 
is  calculated  as  the  magnitude  of  the  right-hand  side  of  eq.  IT  because,  for 
lossy  dielectric  samples  ,  C  will  include  a  contribution  due  to  the  loss  and 
will  provide  a  better  estimate  for  32^Z  than  might  be  obtained  using  eq.  9. 

On  certain  usually  rare  occasions,  the  sample  may  turn  out  to  be  very 
nearly  an  odd  multiple  quarter-guide-wavelength  in  the  sample  material.  In 
this  case,  zQ  is  nearly  a  quarter-guide-wavelength,  so  tan  U  =  tan  (  2tt  zQ/Xg) 
takes  on  very  large  values ,  and  subsequent  calculations  involving  the  hyper- 
bolic functions  of  a2d  and  &2d  exceed  the  capacity  of  the  computer.  There- 
fore, a  test  was  inserted  in  the  program,  and,  if  tan  U  exceeds  a  value  of 
200,  it  branches  to  a  later  point  in  the  program,  where  an  approximate  value 
for  the  dielectric  constant  is  calculated  using  eq.  8  and  the  value  of  32<i 
arrived  at  by  the  FUNCTION  FINDX  subroutine  to  be  described  later.     A  state- 
ment then  appears  in  the  printout,  "Bad  sample  length  for  this  frequency, 
nearly  an  odd  multiple  quarter-wavelength,"  and  zeros  are  printed  out  for 
the  loss  tangent,  loss  factor,  and  conductivity. 

In  the  usual  case  where  tan  U  does  not  exceed  a  value  of  200,  the  esti- 
mate of  32^  is  then  improved  by  adjusting  it  to  satisfy  (tan  32d)/32<i  =  C  in 
the  subroutine,  FUNCTION  FINDX.  This  value  of  B2d  and  a  value  for  a2d,  cal- 
culated from  eq.  15  using  this  improved  estimate  of  &2d,  then  comprise  the 
first  estimate  for  y2^.  These  unadjusted  values  for  a2d  and  $2d  are  printed 
out  as  diagnostics  immediately  following  the  first  estimate  of  32^.  The  com- 
ponents a2d  and  $9d  are  then  adjusted  in  a  subroutine  called  ADJUST,  to  be 
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explained  later,  until  eq.  17  is  satisfied  to  within  a  desired  limit.  The 
adjusted  values  of  a  2d  and  $2d  are  also  printed  out  along  with  the  unadjusted 
values . 

f  sce^len 

Since  (tan  32^)/32^  an(^  (tanh  12^/^2^-  are  J*"'"1  M  ^'r^TTpg  functions  and  each 
solution  to  eqs.  7  and  9  results  in  a  different  set  of  values  for        and  , 
the  program  finds  three  sets  of  a2^Z  and  B2<i  satisfying  eqs.  7  or  17  from 
which  to  proceed  with  the  calculations.     The  relative  permittivity  or  dielec- 
tric constant  and  the  relative  dielectric  loss  factor  of  the  sample  are  then 
calculated  from  an  expansion  of  eq.  6, 
(a2d  +  Q$2d)  A  ]2 

-  jV'  =  P  -  I  ^  -J     =  P  -  S2[(a2d)2  -  (32<i)2]  +  J  S2  2  a2d  $2d, 

where  P  =  (A0AC)2  and  S  =  \Q/2vd.     The  associated  loss  tangent  for  each  of 
the  three  sets  of  values  for  e£  and  e£  is  calculated  (tan  6  =  &^/&').     The  cor- 
rection for  wall  losses  with  the  sample  present  in  the  guide,  tan  6WS ,  is  then 
calculated  using  eq.  13,  after  getting  tan  6^  using  eq.  11,  and  tan  &Ws  is  sub- 
tracted from  tan  6  to  give  the  corrected  value  for  tan  6.    Using  this  corrected 
value,  a  corrected        is  found  as  the  product  of        and  the  corrected  tan  6. 
The  three  sets  of  dielectric  properties  are  printed  out  along  with  important 
input  data  and  a  few  diagnostics,  unless  this  printout  is  suppressed  "by  in- 
structions on  the  parameter  card,  that  is,  nothing  punched  in  columns  68  to 
73.     The  program  selects  the  probable  value  of       ,  that  which  is  closest  to 
the  estimated  value  for        from  the  input  data,  and,  at  the  end  of  the  calcu- 
lations for  a  given  data  set,  prints  out  a  summary  table  including  description, 
date,  moisture  content,  temperature,  sample  identification,  frequency,  and 
values  for      ,  tan  6,  and  conductivity,  calculated  as  a  =  oj  eQ  z£  = 

0.55668/  z'l  for  a  in  millimhos/cm-1 .     Examples  of  the  input  card  data  and 
printout  for  rectangular  and  cylindrical  waveguide  and  coaxial  line  cases  are 
shown  in  figs.  5,  6,  and  7-     The  summary  table  printout  is  shown  in  fig.  8. 

It  should  be  noted  that  the  relationship  of  eq.  2  based  on  twice -minimum- 
power  or  3-dB  points  is  used  in  the  derivation  of  eqs.  10  and  11.  Therefore, 
in  this  program,  the  wall-loss  correction  for  the  sample-filled  part  of  the 
waveguide,  tan  5^s  of  eq.  13,  is  also  based  on  the  assumption  of  a  3-dB  air- 
node  measurement.     Since  the  program  automatically  converts  the  dB  level  for 
the  air-node  measurement  to  that  used  for  the  sample-node  measurement ,  sample- 
node  measurements  on  very  low-loss  materials,  where  the  wall-loss  correction 
is  important,  should  probably  be  performed  using  the  3-dB  level.     In  practice, 
the  need  for  less  than  3-dB  sample-node  measurements  was  encountered  only  with 
high-loss  materials  where  the  wall-loss  correction  was  insignificantly  small 
anyway.     Therefore,  the  3-dB  level  requirement  for  sample-node  measurements 
is  only  necessary  when  working  with  very  low-loss  materials. 
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DIELECTRIC  PROPERTIES  FOR  1970  SCCUT  66  HRW  WHEAT,  X-B AND 
RUN       1  7/30/71 


MOISTURE  CONTENT 

TEMPERATURE    (DEGREES  FAHRENHEIT) 

FREQUENCY  (GIGAHERTZ) 

WAVEGUIDE   HORIZONTAL  DIMENSION  (MM) 

WAVEGUIDE  VERTICAL  DIMENSION  (MM) 

SAMPLE  LENGTH  (MM) 

AIR  NODE  POSITION  (MM) 

SAMPLE  NODE  POSITION  (MM) 

AIR  NODE  WIDTH   (MM)  6.0  CB 

SAMPLE  NODE   WIDTH   (MM)  3.0  Ce 

REFERENCE  DIMENSION  (MM) 
ESTIMATED  DIELECTRIC  CONSTANT 
FIRST  ESTIMATE  OF   B2D  (ESTX) 


10, 
76 

8, 
22, 
10, 
34. 
13, 
22, 

O, 

6, 
130, 

2, 


3880 
8600 
1600 
8000 
6875 
2800 
2910 
3000 
OOOO 
7000 
7042 


UNADJUSTED  A2D,  B2C 
ADJUSTED  A2D,  B2D 


0.33920 
0.42670 


5.73540 
5.82789 


SAMPLE  3-3 


DIELECTRIC  CONSTANT 
LOSS  TANGENT 
LOSS  FACTOR 
CONDUCTIVITY 


RESULTS 


1.512730 
0.087389 
O.  132196 
0.617282 


2.624163 
0.091700 
0.240635 
1. 123630 


4. 204988 
0.075326 
0.316744 
1.479012 


DIELECTRIC   PROPERTIES  FOR     1970  SCOUT  66  HRW  WHEAT,  X-BAND 


DATE 


MOIST  TEfP 
CONT.  (F.) 


SAMPLE 
I  CENT. 


SUMMARY  TABLE 

FREO.  DIELECTRIC  LOSS  LOSS  CONDUCTIVITY 
(GHZ)     CONSTANT     FACTOR  TANGENT  MlLLIMHOS/CM 


7/30/71 


10.5 


76  3-3 


8.388 


2.6242     0.2406  0.0917 


1. 1236 


Figure  5. — Parameter  card  and  sample  data  card  input  with  printout  for  a 
rectangular  waveguide  measurement . 
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DIELECTRIC  PROPERTIES  FOR  1970  SCOUT  66  HRW  WHEAT.  CRL  DIELECTR, 
RUN       1  8/3/71 


MOISTURE  CONTENT 

10. 

5 

TEMPERATURE    (DEGREES  FAHRENHEIT) 

76 

FREQUENCY  (GIGAHERTZ) 

8. 

5144 

WAVEGUIDE  HORIZONTAL  DIMENSION  (MM) 

25. 

4500 

WAVEGUIDE   VERTICAL  DIMENSION  (MM) 

O. 

O 

SAMPLE  LENGTH  (MM) 

46. 

OOOO 

AIR  NODE  POSITION  (MM) 

61. 

1675 

SAMPLE  NODE  POSITION  (MM) 

79. 

9600 

AIR  NODE  WIDTH   (MM)                  3.0  DB 

O. 

0890 

SAMPLE  NODE  WIDTH   (MM)           3.0  CB 

9. 

3200 

REFERENCE  DIMENSION  (MM) 

300. 

OOOO 

ESTIMATED  DIELECTRIC  CONSTANT 

2. 

7000 

FIRST  ESTIMATE  OF   E2D  (ESTX) 

8. 

5951 

UNADJUSTED  A2C,  B2C 
ADJUSTED  A2D,  B2D 


0.42110 
0.61110 


8.67469 
8.69968 


SAMPLE  3-1 


DIELECTRIC  CONSTANT 
LOSS  TANGENT 
LOSS  FACTOR 
CONDUCTIVITY 


RESULTS 


1.773952 
0.088833 
0.157586 
0.746929 


2.641019 
0.077709 
0.205232 
0.972763 


3.794526 
O. C59573 
O. 226052 
1.071445 


DIELECTRIC  PROPERTIES  FOR     1970  SCOUT  66  HRW  WHEAT,  CRL  DIELECTR. 


DATE 


8/3/71 


MOIST  TEMP 
CONT.  (F.) 


SAMPLE 
I  CENT. 


10.  5 


76  3-1 


SUMMARY  TABLE 

FREQ.  DIELECTRIC  LOSS  LOSS  CONDUCTIVITY 
(GHZ)     CONSTANT     FACTOR  TANGENT  MILLIMHOS/CM 


8.514       2.6410     0.2052  0.0777 


0.9728 


Figure  6. — Parameter  card  and  sample  data  card  input  with  printout  for  a 
cylindrical  waveguide  measurement,  using  guide-wavelength  data  for  the 
frequency  input. 
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DIELECTRIC  PROPERTIES  FOR  1970  SCOUT  66  HRW  WHEAT,  RES  SLRC-LMC 
RUN       1  8/2/71 


MOISTURE  CONTENT 

TEMPERATURE   (CEGREES  FAHRENHEIT) 

FREQUENCY  (GIGAHERTZ) 

WAVEGUIDE  HORIZONTAL  DIMENSION  (MM) 

WAVEGUIDE  VERTICAL  DIMENSION  (MM) 

SAMPLE  LENGTH  (MM) 

AIR  NODE  POSITION  (MM) 

SAMPLE  NODE  POSITION  (MM) 

AIR  NODE  WIDTH  (MM)  3.0  CB 

SAMPLE  NODE  WIDTH   (MM)  3.0  OB 

REFERENCE  DIMENSION  (MM) 
ESTIMATED  DIELECTRIC  CONSTANT 
FIRST  ESTIMATE  OF   E2D  (ESTX) 


UNADJUSTED  A2D,  B2D 
ADJUSTED  A2D,  B2D 


O. 16494 
O. 15744 


10.5 
76 

2.4356 

0.0 

0.0 
38.3000 
94.7500 
71.0100 

O.IOOO 

3. 8200 
150.0000 

3.0000 

3.3873 


3.30719 
3.17719 


SAMPLE  3-2 


RESULTS 


DIELECTRIC  CONSTANT 
LOSS  TANGENT 
LOSS  FACTOR 
CONDUCTIVITY 


2.632839 
0.098944 
0.260504 
0.353200 


10.549779 
0.101583 
1.071676 
1.453012 


23.807678 
0.105471 
2. 511011 
3.404506 


DIELECTRIC   PROPERTIES  FOR     1970   SCOUT  66  HRW  WHEAT,    R£S  SLRC-LMC 


DATE 


SUMMARY  TABLE 

MOIST  TEMP  SAMPLE  FREQ.  DIELECTRIC  LOSS  LOSS  CONDUCTIVITY 
CONT.    (F.)        IDENT.        (GHZ )     CONSTANT     FACTOR  TANGENT  MILLIMHOS/CM 


8/3/71  10.5 


76  3-2 


2.436       2.6328     0o2605  0.0989 


0.3532 


Figure  7. — Parameter  card  and  sample  data  card  input  with  printout  for  a 
coaxial  line  measurement,  using  guide-wavelength  data  for  the  frequency 
input . 
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DIELECTRIC   PROPERTIES  FOR     1970  SCOUT  66  HRW  WHEAT,  X-BAND 


SUMMARY  TABLE 


DATE 


MOIST  TEMP 
CONT.  (F.) 


7/30/71  10.5 

7/30/71  10.5 

7/30/71  10.5 

7/30/71  10.5 


SAMPLE 
I  CENT. 


76  3-3 

76  3-3 

76  3-3 

76  3-3 


FREO.  DIELECTRIC  LOSS  LOSS  CONDUCT  I VI TY 
(GHZ)     CONSTANT     FACTOR  TANGENT  MILLIMHOS/CM 


8.388  2.6242 

10.186  2.6005 

11.030  2.5864 

12.228  2.5174 


0.2406  0.0917  1.1236 

0.2582  0.0993  1.4643 

0.2784  0.1076  1.7095 

0.2143  0.0851  1.4590 


Figure  8. — Parameter  card  and  sample  data  card  input  with  summary  printout 
only  for  series  of  rectangular  waveguide  measurements. 


Program  Subroutines 


FUNCTION  FINDX 


In  the  main  program,  an  estimate  for  is  obtained  using  eq.  8.  The 

FORTRAN  FUNCTION  FINDX  is  used  to  solve  (tan  62^)/(32^)  =  C  for  when  C 

is  known  from  measurements  data.     In  the  program,  X  is  used  to  represent  32^- 
Because  of  the  nature  of  the  (tan  X)/X  function,  some  care  must  be  exercised 
to  obtain  desired  solutions  by  iterative  techniques.     A  plot  of  this  function 
is  shown  in  fig.  9  for  the  range  0  £  X  ^  Stt  .     For  (tan  X)/X  =  C,  it  is  ob- 
vious that  no  solutions  can  exist  in  the  first  quadrant  for  C  <  1.     For  a 
given  value  of  C,  there  are  many  values  of  X  which  will  satisfy  the  equation. 
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Figure  9- — Graph  of  the  function  (tan  X)/X. 

The  program  calculates  and  prints  out  dielectric  properties  correspond- 
ing to  three  different  values  for  X,  or  62^5  one  nearest  the  value  derived 
from  the  estimated  dielectric  constant  of  the  input  data,  and  one  on  either 
side  of  this  value  if  possible.     If  the  estimated  value  of  X  obtained  from 
eq.  8  is  greater  than  it,  tt  is  subtracted  from  this  value  to  obtain  the  esti- 
mated X  value  on  the  lower  side  of  the  original  X  estimate.     If  this  results 
in  a  value  for  X  in  the  first  quadrant,  the  program  checks  to  see  whether 
C  ^  1,  and,  if  it  is,  tt  is  added  back  on  to  take  the  estimated  X  into  the 
third  quadrant  where  a  possible  solution  exists.     If  the  estimated  X  does 
not  lie  in  the  first  quadrant,  or,  if  it  does,  but  C  >  1,  this  value  is  taken 
for  the  first  set  of  calculations  and  stored  for  printout  as  the  estimated 
diagnostic  immediately  following  the  printout  of  input  data.     The  pro- 
gram then  sets  32^  =  FINDX  which  is  a  function  of  C  with  proper  sign  and  the 
estimated  value  of  X. 


The  FUNCTION  FINDX  subroutine  finds  the  value  of  X  satisfying 
(tan  X)/X  =  C  in  the  following  way:     First,  the  quadrant  containing  the  esti- 
mated X  is  identified  and  appropriate  limits  are  established  within  which  the 
search  for  the  proper  value  of  X  will  be  carried  out.     A  test  is  also  made  to 
determine  whether  C  and  (tan  X)/X  have  the  same  sign.     If  they  do  not,  and 
the  value  of  X  lies  near  an  asymptote  (ascertained  by  |c|  >  2),  the  estimated 
X  is  shifted  to  the  other  side  of  the  asymptote  a  distance  equal  to  its  dis- 
tance from  the  asymptote  (points  a  and  b,  fig.  9)-     This  will  change  the  sign 
of  X  thereby  giving  C  and  (tan  X)/X  the  same  sign  so  that  a  solution  will  be 
possible.     New  limits  for  the  search  are  also  then  established  (tt/2  and  Sir / 2 
for  point  b,  for  example).     If  the  sign  of  C  and  that  of  (tan  X)/X  do  not 
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agree  and  the  estimated  X  is  not  near  an  asymptote  (|c|  ^  2),  the  estimated 
X  is  shifted  to  the  other  side  of  nu,  in  which  quadrant  a  solution  can  "be 
found  (points  o  and  d,  for  example). 

An  iterative  technique  employing  Newton's  successive  approximation 
method,  eq.   l6 ,  is  then  used  to  find  the  proper  value  of  X  satisfying 
(tan  X)/X  =  C.     Here  f(x)  =  (tan  X  -  CX)/X  and  f (x)  =  (X  sec2  X  -  tan  X)/X2 
so  that  f{x)/f'(x)  =  Z(tan  X  -  CX)/{X  sec2  X  -  tan  X) .     It  may  happen  that 
the  value  of  J2  obtained  using  the  estimated  X  value  for  X\  in  eq.  l6  will 
fall  outside  the  prescribed  limits.     For  example,  a  solution  is  needed  "be- 
tween A  =  it/2  and  B  =  3n  /  2  in  fig.  10  and  the  first  value  obtained  for 


4  8 
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*, 

Xa 

r 

I  7^  *  

■  y 

r 

Figure  10. — Graph  of  the  function  (tan  X)/X  -  C  showing  limits 
for  successive  approximation  solution. 

X2  >  3tt / 2 .     The  program,  therefore,  rejects  X^  outside  the  prescribed  region 
and  sets  the  new  estimate  of  X  equal  to  (X    +  Sir/ 2)/ 2  and  proceeds  until  X2 
falls  within  the  prescribed  limits.     Successive  approximations  then  proceed 
until  the  difference  between  two  successive  approximations  of  X  is  less  than 
0.0001  and  the  program  then  sets  the  last  value  for  X  =  FINDX. 

ADJUST 

Following  the  use  of  the  FUNCTION  FINDX  subroutine  to  improve  the  esti- 
mate of  32<i  and  the  subsequent  calculation  of  an  estimate  for  a„J,  the  ADJUST 
subroutine  is  used  to  adjust  the  values  of  the  components  of  y =  a^d  +  jQ^d 
to  satisfy  eq.   17.     Here  UQ  and  VQ  represent  the  real  and  imaginary  parts  of 
the  right-hand  side  of  eq.  IT,  and  their  values  are  determined  by  measurement 
data.     The  subroutine  AC0MP  provides  values  for  (tanh  y2d)/y2d  =  Va  +  Va 
from  the  initial  a2d  and  $2d  values  calculated.     The  ADJUST  subroutine 
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I  START  1 
I 

|delta=  o.oi  I 


CALCULATE   UA, VA 


|  DIFF  I  | 


CALCULATE  UA, VA  I 

  1  

|  PI  FT  2 


|  DELTA  ^DELTA/2  [ 


Figure  11. — Flow  chart  for  the  ADJUST  subroutine, 


-  21  - 


alternately  adjusts  values  of  a^d  and  Q^d  by  an  increment,  delta,  and  tests 
to  see  whether  the  resulting  V a  +  jVa  is  approaching  the  desired  value  of 
U    +  jV  .     This  is  done  by  calculating  the  difference  in  the  two-space  vec- 
tors, U    +  jVQ  and  Ua  +  jVa,  /(f/0  -  Ua)2  +  (VQ  -  Va)2,  and  requiring  that 
this  value  decrease  as  adjustments  are  made.     When  no  further  improvement 
can  be  found  by  either  positive  or  negative  additions  of  delta  to  a2d  and 
32<^5  the  value  of  delta  is  halved,  and  the  adjusting  procedure  repeated. 
This  reduction  of  delta  and  adjusting  of  a2d  and  &2d  values  continues  in 
accordance  with  the  flow  chart  of  fig.  11  until  the  magnitude  of  the  vector 
difference  is  less  than  0.0002,  or  until  the  delta  value  has  been  halved  12 
times,  when  it  returns  to  the  main  program.     If  the  vector  difference  is  not 
reduced  to  a  value  less  than  0.0002,  the  value  of  the  difference,  called 
DIFF1  in  the  program,  is  printed  out  and  identified  with  one  of  the  three 
dielectric-constant  determinations  to  which  it  applies. 
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APPENDIX 


Computer  Program  Listing 


FORTRAN    IV   G  LEVEL      19  MAIN  DATE  =  72004  14/28/14 

C  LATEST  REVISION  MADE  AUGUST   10,    1971.  SHORTLIN 

C  MAIN  PROGRAM  FOR  CALCULATING  DIELECTRIC   PROPERTIES  FROM 

C  SHORT-CIRCUITED  COAXIAL-LINE  OR   WAVEGUIDE  MEASUREMENTS 

C  INPUT   r ATA  DESCRIPTION 

C  LL                  NUMBER  OF   SETS  OF  DATA 

C  JJ                  NUMBER   OF  READINGS   PER  SET 

C  SAMPL            DESCRIPTION  OF  THE  DATA  SET,     40  SPACES 

C  HDIMW           HORIZONTAL  DIMENSION  OF   THE   WAVEGUIDE    IF  RECTANGULAR 

C  DIAMETER  OF  THE   WAVEGUIDE    IF  CIRCULAR 

C  ZERO  FOR  COAXIAL  WAVEGUIDE 

C  VDIMW           VERTICAL   DIMENSION  OF  THE    WAVEGUIDE    IF  RECTANGULAR 

C  ZERO  OTHERWISE 

C  RDIMW           WAVEGUIDE  REFERENCE   DIMENSION,    MEASURED   FROM  THE 

C  SHORT  TO  THE  0  MM.    REFERENCE   POINT  OF   THE  PROBE. 

C  SLOT              WIDTH  OF   SLOT   IN  SLOTTED  SECTION 

C  ALCALC         ANY  POSITIVE  NUMBER   IF  COMPLETE   RESULTS   ARE  TO  BE 

C  PRINTED 

C  ZERO   IF   SUMMARY  TABLE  ONLY   IS  DESIRED 

C  SAMID            INDIVIDUAL   SAMPLE  DESCRIPTION,    8  SPACES 

C  DATE             E.G.  12-25-70 

C  MC                 MOISTURE  CONTENT  OF  SAMPLE,  4  SPACES 

C  TEMPF           SAMPLE  TEMPERATURE    IN  DEG.    FAHRENHEIT,   4  SPACES 

C  FREQ             FREQUENCY   IN  GIGAHERTZ 

C  ANODP           AIR  NODE   PLUS  HALF  NODE  WIDTH 

C  ANODM           AIR  NODE   MINUS   HALF  NODE  WIDTH 

C  SNOOP           SAMPLE  NODE   PLUS   HALF  NODE  WIDTH 

C  SNODM           SAMPLE  NODE  MINUS  HALF  NODE  WIDTH 

C  DDI  MS           DEPTH  OF  THE  SAMPLE 

C  ADB               OB   LEVEL   FOR  AIR-NODE  MEASUREMENT 

C  SDB                DB  LEVEL   FOR   SAMPLE-NODE  MEASUREMENT 

C  ESTDC           ESTIMATED  DIELECTRIC  CONSTANT 

C  WL                  WAVELENGTH   IN  FREE  SPACE 

C  WLC               CUTOFF  WAVELENGTH 

C  WLG                WAVELENGTH   IN  THE  WAVEGUIDE 

C  ALL   MEASUREMENTS  ARE   IN  MILLIMETERS 

C 

0001  REAL  MC,LTW,LTWS 

000  2  DIMENSION  ALPHA (3 ) , BETA ( 3 ) , DC ( 3 ) , ALT ( 3 ) .FACTOR (3 )  , SDC (50)  ,  SLT ( 50  )  , 
*  SFACT(50),SCOND(50) , SF ( 50 ) , ST F ( 50 ) , S AMPL ( 10 ) , S AM  ID ( 2  )  ,  0ATE(2», 
*SDAY(50,2),SSMP(50,2>,SMC(50),DIFF(3),COND(3J 

0003  I  R  E  AD=  5 

0004  IRITE=6 

0005  TNU=0 

0006  ESTEPM=2.5 
C 

0007  PI=2. 141592654 

0008  REAC(IREAC,300)  LL 

0009  300  F0PfAT(I5> 

0010  DO   500  L=1,LL 

0011  WRITE(IRITE,1200) 

0012  1200  FORMAT(lHl) 

0013  REAC ( IREAD , 305 )    JJ,SAMPL,HDIMW, VDIMW, RDIMW, SLCT, ALCALC 

0014  305   FORMAT! I  3  , 10 A4 , 5 F6. 2 ) 

0015  SL0TC=1 

0016  IF    (SLOT. ECO. AND. HDIMW. GT.O. AND. VDIMW. GT.O)    SLOT  =  3.20 

0017  IF    (SLOT. ECO. AND. HDIMW. GT.O.  AND. VDIMW.  EQ.C)  SL0T=1.24 

0018  IF    (HDIMW. GT.O. AND. VDIMW. GT.O)  WLC=HDIMW*2 

0019  IF    (HDIMW. GT.O. AND. VDIMW. EO.O)    WLC=HDIMW*  1.706293 
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0020  DO  400  J=1,JJ 

0021  DAT A=0 

00  22  22  REACdREAD, 311  »SAMID,  DATE, MCTEMPF,  FRECANODP,  ANOQM,  SNGDP,  SNODM,  DD 
*IMS  ,ADB, SCB, ESTDC 

0023  311  F0f^ATt2A4,2A4,A4,lX,A4,F7.4,5F7.3,2F4.0,F3.11 

0024  TF=TEMPF 

0025  IF    (FREQ.GT.O.AND.DDIMS.GT.O)   GO  TO  56 

0026  IF    (FREQ.EQ.O.AND.DDIMS.GT.O)   GO  TO  54 

C 

0027  50   WLG=ABS ( ANC DP+ANO DM-SNODP-SNODM I 

0028  IF    (HDIMW.EQ.O)    GO  TO  53 

0029  IF   (VDIMW.EQ.O)   GO  TO  43 
C030  WLGN=WLG 

0031  51  CUBIC=SL0T**2*WLGN**3«-8*PI*VDIMW*HDIMW**3*(WLGN-WLG) 

00  32  CUBPM=3*SL0T**2*WLGN**2+8*PI*VDIMW*HDI MW**3 

0033  GO  TO  45 

0034  43  WLGN=WLG 

00  35  FOFAW= 1.41 842*SL0T**2/( 8*1 PI*HDIMW/2 )**4) 

0036  44  CUBIC*F0FAW*WLGN**3+WLGN-WLG 

0037  CUBPM=3*F0FAW*WLGN**2+1 

0038  45  WLGN1=WLGN-CUBIC/CUBPM 

0039  IF    (ABS(WLGNl-WLGN).LE. 0.0001)    GO  TO  52 

0040  WLGN=WLGN1 

0041  IF    (VDIMW.EQ.O)   GO  TO  44 

0042  GO  TO  51 

0043  52   WLG  =  WLGM 

0044  53    IF    (HDIMW.GT.O)    WLG= WLG*WLC / SORT ( WLG**2 *W LC**2 ) 

0045  F  =  2<59. 6S66/WLG 

0046  GO  TO  22 

0047  54  FREQ=F 

0048  56  F=FREQ 

0049  IF    (FREQ.GT. 0.1. AND.FREQ.LT. 13)    GO  TO  58 

0050  315  DATA=1 

0051  GO  TO  350 

0052  58   WL  =  299.6966/FRE0 

0053  IF    (HDIMW.EQ.O. AND. VDIMW.EQ.O)  WLG=WL 

0054  IF    (HDIMW.GT.O)    WLG=W L /SORT ( 1- ( WL / WL C ) **2 I 

0055  IF    (HDIMW.EQ.O)   GO  TO  86 

0056  IF    (VDIMW.EQ.O)    GO  TO  59 

0057  SL0TC=1*( SL0T*WLG)**2/(8*PI*VDIMW*HDI MW**3 ) 

0058  GO  TO  66 

0059  59  C0N=1. 84118 

006  0  SL0TC=1*( SLOT*WLG*CON )**!/ ( 8*( PI*HDIMW/2) **4* (C0N**2»1 )  ) 

0061  GO  TO  86 

0062  86  PRNT1=ABS( ANODP-ANODM) 

0063  PRNT2  =  ABS( SNODP-SNCDM  ) 

0064  ZANODE=( ANOO P+ ANODM ) / 2 

0065  ZSNCDE=( SNGDP+SNODM ) /2 

0066  DLTZA=ABS(ANODP-ANODN)/SLOTC 

0067  DLTZS=ABS( SNOOP -SNODM) /SLOTC 

0068  87  SHIFT=(ZAN0DE-ZSNODE)/SL0TC 

0069  DE=ZANODE*RDIMW 

0070  DET=0 

0071  88  DET=DET+WLG/2 

0072  IF    (ABS ( DE-DET ) .GE.WLG/4)   GO  TO  88 

0073  DE=DET 

0074  IF    (ACB.EC.O)  ADB=3 

0075  ARG=0.2302585*ADB 
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0076 
0077 
0078 
0079 
0080 
0081 
0082 

ooe3 

0084 
0085 
0086 
0C87 
0088 
0089 

0090 
0091 
0092 
0093 
0094 
0095 
0096 
0097 
0098 
0099 
0100 
0101 
0102 
0103 


0104 
0105 


0106 
0107 
0108 
0109 
0110 

0111 
0112 


0113 
0114 
0115 
0116 
0117 
0118 
0119 
0120 
0121 
0122 
0123 
0124 
0125 


P0WER=EXP(ARG) 

VSWR=SQRT( P0WER-C0S< P I *DLTZA/ WL G ) **2 )/SIN(PI*DLTZA/WLG) 
IF    (SDB.EQ.01  SDB=3 
ARG*0.2302585*SDB 
POW  ER=E  X  P ( ARG ) 

ARG=SQRT ( (POWER-1  ) / ( VSWR**2- 1 ) I 
DLTZA=WLG*ARSIN(ARG)/PI 
RAT  10= ( DE-DOI MS-SHI  FT) /OE 
95  PHIPR=PI*(DLTZS-RATIO*DLTZA)/WLG 

VSWP  =  SORT ( POWER-COS ( PHI  PR) **2 ) /SIN( PHI  PR) 
PHI=1/VSWP 

100  ZO=-(SHIFT+DDIMS) 

IF    (ZO.GE.O)   GO  TO  102 

101  Z0=Z0*WLG/2 

I.E.      ADD  HALF-WAVELENGTHS  UNTIL   ZO  BECOMES  POSITIVE. 
IF    (ZO)  101,102,102 

102  U=2*PI*ZC/WLG 
TNU=TAN ( U ) 

IF    (HDIMW.GT.O)    P = ( WL /WLC ) ** 2 

IF    (HDIMW.EQ.O)  P=0 

Q=WLG/(2*FI*QDIMS ) 

S  =   WL/( 2*PI*DDIMS) 

DEN=  0/(l*(PHI*tNUJ**2) 

A=TNU*(PHI**2-1 ) 

B=PHI*( 1+TNU**2) 

C  =  SORT (A*A*B*B)*DEN 

C1=C*SIGN< 1.0, A) 

UO=A*D€N 

VO=-B*DEN 

ESTX  IS  THE  FIRST  ESTIMATE  OF  B2D,   COMPUTED  FROM  E •    FOR  THE 
APPROPRIATE  GUIDE  WITH  A2D=0. 
IF    (ESTDC.EQ.O)   EST  DC=  ESTEPM 
ESTX  =  SORT(ESTDC-P)/S 

KWAO  CORRESPONDS  TO  THE  QUAORANT  INTO  WHICH  ESTX   (IN  RADIANS)  FALLS, 

104  KWAD=  1 

IF    (ESTX. GT. PI)   ESTX  =  ESTX-PI 
2  IF    (ESTX.LT.KWAD*PI/2)   GO  TO  5 
KWAC=KWAD*1 
GO  TO  2 

SOLUTICNS  DO  NOT  EXIST   IN  THE  FIRST   QUADRANT    IF  C.LT.l. 

5  IF    (KWAD-1)  6,6,40 

6  IF    (Cl.LE.l)  ESTX=ESTX*PI 

FIND  THREE  VALUES  OF  B2D  AND  A2D. 
40  ESTXO=ESTX 

DO  42  1=1,3 
42  DIFF(I)=0 

IF    (ABS(TNU) .GT.200)   GO  TO  230 

DO   110  1=1,3 

B2D=FINDX(C1, ESTX ) 

ESTX=B2D*PI 

105  A2D=(Q*B2D*PHI*( 1+TNU**2) )/( ( 1+TAN ( B2D ) **2 J-TAN ( B2 D ) /B2C ) 
IF    (1-1)  106,106,1^7 

106  TB2C=B2D 
TA2C=A2D 

107  CALL  ADJUST   ( A2 D, B2 D, UO, VO , I , D I F ) 
IF    (DIF.LT. 0.0002 )  GO  TO  108 
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0126 

DI  FF ( I ) =DI F 

0127 

108 

ALPHA ( I )=A2D 

0128 

BETA( I )=B2D 

0129 

110 

CONTINUE 

01  30 

c 

116 

IF    (HDIMW.EQ.O)   GO  TO  200 

c 

CALCULATE  3  DIELECTRIC  CONSTANTS,   LOSS  FACTORS,  AND  LOSS 

c 

FOR  CIRCULAR  OR  RECTANGULAR  WAVEGUIDES, 

0131 

DO    140  1=1,3 

0132 

DC(I)   =   P-S**2*(ALPHA( I)**2-BETA(I  )**2> 

0133 

FACTOR(I)   =  S**2*2*ALPHA( I ) *BETA( I  J 

0134 

ALT(I)=FACTOR(I  )/DC(I) 

0135 

LTW=(1-P)*CLTZA/DE 

0136 

IF    (VDIMW.EQ.O)    GO  TO  132 

0137 

LTWS=LTW*( HDI MW/( 2*VDIMW)*P/DC (I ) )/(HDIMW/(2*VDIMW)+P) 

0138 

GO  TO  135 

0139 

132 

LTWS=LTW*(0.42+P/DC (I )  )/(0.42«-P) 

0140 

135 

ALT ( I  )  =ALT ( I )-LTWS 

0141 

FACTOR ( I )=ALT(I)*DC(I ) 

0142 

140 

COND(I)   =  .55668*FACTOR(  H*FREO 

0143 

GO  TO  350 

c 

C     CALCULATE  3  DIELECTRIC  CONSTANTS,    LOSS  FACTORS,   AND  LOSS  TANGENTS  FOR 


C 

COAXIAL  LINES. 

0144 

200 

DO  220  1=1,3 

0145 

DC( I)=S**2*(BETA( I )**2-ALPHA( I )**2) 

0146 

FACTOR ( I  )  =  S**2*2* ALPHA! I ) *BETA ( I ) 

0147 

ALT( I )=FACTOR(I »/DC(I ) 

0148 

LTWS=DLTZA/DE 

0149 

ALT(I)=ALT(I)-LTWS 

0150 

FAC70RU  )=ALT(  I)*DC(I  ) 

0151 

220 

COND(I)    =  .55668*FACT0R(I )*FREQ 

0152 

GO  TO  350 

0153 

230 

DATA=1 

0154 

DO   231  1=1,3 

0155 

B2D=FINDX(C1,ESTX) 

0156 

DC(  I)=P+(WL*B2D/(2*PI*DDIMS> )  **2 

0157 

231 

ESTX=B2D+PI 

0158 

DO  232  1=1,3 

0159 

ALT ( I ) =0 

0160 

232 

FACTOR ( I )=0 

C 

C     PRINT  COMPLETE  RESULTS  IF  ALCALC  IS  OTHER  THAN  ZERO. 

0161  350  IF   (ALCALC. EQ.O)   GO  TO  380 

0162  WRITE(IRITE,375)   SAMPL , J , DAT E, MC 

0163  375  FOR  MAT( 1H0 , •      DIELECTRIC   PROPERTIES   FOP   •  ,10A4,//3X, 'RUN' ,I4,5X,2A 

*4///,*   MOISTURE  CONTENT*  »  31X, A4) 

0164  WRITE( IRITE.376)   TF , F , HDI MW , VDI MW , DDI  MS , Z ANO DE , ZSNODE  ,  ACB.PRNTl, 

*  SDB, PRNT2.RDI MW, ESTDC 

0165  376  FOR  MAT  (  1H   ,  'TEMPERATURE   (DEGREES  FAHRENHEIT)-'  ,12X,A5/'   FREQUENCY  ( 

♦GIGAHERTZ) • ,F33. 4, /•    WAVEGUIDE   HORIZONTAL   DIMENSION    ( MM ) ' , Fl 9.4,  / • 

*  WAVEGUIDE  VERTICAL  DIMENSION   (MM)  •  ,F21.4, /'    SAMPLE   LENGTH  (MM)*,F 
♦36.4, /•    AIR  NODE  POSITION   (MM)         «,F28.4,/*   SAMPLE  NODE  POSITION 
*(MM) • ,F29.4/*   AIR  NODE  WIDTH   (MM)»,F11.1,      •    De     *,F19.4,      /■  SAMP 
*LE  NODE   WIDTH   (MM)*,F   8.1,'   DB*,F21.4,/*  REFERENCE  DI  MEN 
*SION  (MM)« ,  F2-0.4, /■    ESTIMATED  DIELECTRIC  CONSTANT* , F25.4) 

0166  IF    (DATA.EQ.O.O)   GO  TO  373 

0167  SDC(J)=0 
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0168 
0169 
0170 
0171 
0172 
0173 

0174 
0175 
0176 
0177 
0178 
0179 
0180 
0181 
0182 
0183 
0184 
0185 

0186 

0187 
0188 
0189 
0190 
0191 
0192 

0193 
0194 
0195 
0196 
019T 
0198 
0199 
0200 
0201 
0202 
0203 
0204 
0205 
0206 
0207 
0208 
0209 
0210 
0211 

0212 
0213 
0214 
0215 


0216 
0217 


372 

373 
374 


382 
383 

377 


1100 
379 


SLT(J)=0 
SFACT( J)=0 
SC0ND( J)=0 
EST  EPM=ESTDC 
WRITE(IRITE,371)  TNU 
371   FORMAT ( lHOf ' BAD  SAMPLE   LENGTH  FOR  THIS   FREQUENCY,    NEARLY  AN  OOD  MU 
*LTI PLE  QUARTER  WAVELENGTH. ',/, •   TANGENT ( ( 2*P I *Z0 I / WL G 1  =',F10.3> 
WRITE( IRITE,372)    ( DC ( I >  ,  I  =  1 , 3 ) 

FORMAT ( 1 H   , 'DIELECTRIC  CONSTANT 3F 1 5. 6 , 1H 1 ) 

GO  TO  380 

WRITE( IRITE, 374)  ESTXO 

FORMAT ( 1 H   , 'FIRST   ESTIMATE  OF   E2C   (  E ST  X  )  •  ,  F  26.  4/ ) 
WRITE( IRITE, 382)   TA2D,  TB2D 
WRITEi IRITE, 383)   ALPHA(l),    BETA ( 1 ) 
FORMAT ( 1H   , "UNADJUSTED  A2D,    B2D  '.2F20.5) 
FORMAT ( 1H    « ' ADJUSTED  A2D,   B2D  ',2F20.5//» 
WRITEIIRITE.377ISAMID 

F0RMAT(18X,'      SAMPLE      ',2A4,'      RESULTS  '() 
WRITE! IP ITE  ,378)    (DC (I), 1  =  1, 3), (ALT ( I) ,1  =  1 .3) , ( FACTOR ( I  1,1  =  1,3),  ( 
*CONC( I  )  ,1  =  1,3) 

378  FORMATC    CIELECTRIC   CONSTANT* ,5X,3F15.6,/ '   LOSS  T ANGENT  • , 1 2X ,  3F1 
*5.6,/'    LOSS  FACTOR   « , 1 2X , 3F 1 5. 6 , / •   CONDUCTIVITY* ,12X ,3F15.6) 
DO  379  1=1,3 

IF  (DIFF( I J.EQ.O)  GO  TO  379 
WP ITE( IR ITE, 1100)  DIFF(I),I 

F0RMAT(1H0,'DIFF1=',F10.6,'    FOR   ESTIMATE  NO.', 13/) 
CONTINUE 
WRITE(IRITE,1200) 
ADD  TO  SUMMARY  TABLE   PROBABLE   DC  AND  LOSS-TAN  VALUES 
380  K0UNTR=1 

DELTA=ABS<  CC( 1  )-ESTDC ) 
DO  320  1=2,3 
GUESS=ABS( DC( I ) -EST DC ) 
IF    (GUESS-DELTA)  360,360,320 
KOUIVTR=I 
DEL  TA=GUESS 
CONTINUE 

ESTEPM=  DC ( KOUNTR ) 
SDC(J)=DC(KOUNTR) 
SLT(J)    =  ALT ( KCUNTR ) 
SFACT(J)   =  FACTOR (KOUNTR) 
SCOND(J)    =  ,55668*FACT0R (KOUNTR )*F 
SF(J)    =  F 
SMC( J)=MC 
STF(J)   =  TF 
DO  400  N=l,2 
SDAY(J.N)  =DATE(N) 
SSMP(J.N)    =   SAM  I D (N ) 
WRITE   SUMMARY  TABLE 
405  WRITE(IRITE,410)  SAMPL 

410  FOR MAT( 1H0,' DI ELECTRIC   PROPERTIES   FOR  «,10A4///) 

WRITE( IRITE  ,412) 
412    FORMAT ( 34X, • SUMMARY  TABLE'//'        DATE         MOIST  TEMP 
*.    DIELECTRIC     LOSS       LOSS       CONDUCTIVITY' /11X, 'CO NT. 
*       (GHZ)     CONSTANT     FACTOR  TANGENT  MI L LI MHCS/CM' / ) 
DO  414  J=1,JJ 

414  WRITE (IRITE, 415)    (SDAY(J,N) ,N=1,2),SMC(J) ,STF(J), (SSMP(J,N),N=1,2) 
*,SF(J),SDC(J),SFACT( J),SLT( J),SCOND( J) 


360 


320 


321 


322 


400 


SAMPLE 
(F.) 


FREQ 
I  DENT. 
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0218  415  FORMAT ( 1H   , 2A4, 3< 2X , A4 ) , A4, F7. 3 ,F9.4 , 2F8. 4 ,F 12.4) 

0219  500  CONTINUE 

0220  CALL  EXIT 

0221  END 


FORTRAN  IV  G  LEVEL     19  ADJUST 


0001 

SUBROUTINE  AD JUST ( A2D , B2D, UO , VO , I , DI FF1 J 

0002 

J=0 

0003 

DELTA=0.01 

0004 

CALL  AC0MP(UA,VA,A2D,B2D) 

0005 

DI  FF1  =  SQRT ( (UO-UA 1**2 ♦ ( VO-VA ) **2 ) 

0006 

5 

IF    (DIFFl.LT. 0.0002)  RETURN 

0007 

10 

K=0 

0008 

TA2C=A2D*0ELTA 

0009 

15 

CALL  AC0MP(UA,VA,TA2D,B2D) 

0010 

DIFF2=SQRT( (UO-UA )**2*( VO-VA )**2) 

0011 

K=K+1 

0012 

IF   ( DI FF  2-D IFF  1 )  20,25,25 

0013 

20 

A2D=TA2D 

0014 

22 

DIFF1=DIFF2 

0015 

GO  TO  5 

0016 

25 

IF   (K-l)  30,30,40 

0017 

30 

TA2C=A2D- DELTA 

0018 

GO  TO  15 

0019 

40 

K=0 

0020 

TB2D=82D+DELTA 

0021 

45 

CALL  AC0MP(UA,VA, A2D.TB2D) 

0022 

DI FF2=SQRT  C ( UO-UA )**2  +  ( VO-VA ) **2 ) 

0023 

K=K*1 

0024 

IF   ( DI FF2-DIFF 1 )  50,55,55 

0025 

50 

B2D=TB20 

0026 

GO  TO  22 

0027 

55 

IF   (K-l)  60,60,80 

0028 

60 

TB2  D=B20-DELTA 

0029 

GO  TO  45 

0030 

80 

DELTA=  DELTA/ 2 

0031 

J  =  J+1 

0032 

IF   (J-13)  10,90,10 

0033 

90 

RETURN 

0034 

END 

FORTRAN  IV  G  LEVEL     19  ACOMP 


0001  SUBROUTINE  ACOMP   (UA, VA , A2D, B2D I 

0002  X=A2D*2 

0003  Y=B2D*2 

0004  V1  =  A2D*SINH(X»-*B2D*SIN(Y) 

0005  V2=A2D*SIN(Y)-B2D*SINH(X» 

0006  BOT=( A2D**2+B2D**2r*(C0SH(X)+C0S(Y) ) 

0007  UA=V1/B0T 

0008  VA=V2/B0T 

0009  RETURN 

0010  END 
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FORTRAN  IV  G  LEVEL     19  FINDX 


OOOl 

FUNCTION  FINDX  (CESTX) 

0002 

K=l 

0003 

PI=3»141592654 

0004 

XO=ESTX 

0005 

NSTOP=0 

0006 

10 

IF   ( ESTX.LT.K*PI/2 )  GO  TO  20 

0007 

K=K*1 

0008 

GO  TO  10 

0009 

20 

TEST=TAN(ESTX)/ESTX 

0010 

IF   (C*TEST,GT.0.0 )  GO  TO  21 

0011 

IF   (ABS(C).GT.2.0)  GO  TO  30 

0012 

21 

DO  26  N=l,10 

0013 

IF   (K-(2*N-D)  25,22,26 

0014 

22 

IF    (K-IJ  23,23,24 

0015 

23 

A=0 

0016 

B=PI/2 

0017 

GO  TO  40 

0018 

24 

B=K*PI/2 

0019 

A=B-PI 

0020 

GO  TO  40 

0021 

25 

A=( K-l J*PI/2 

0022 

B=A*PI 

0023 

GO  TO  40 

00  24 

26 

CONTINUE 

0025 

30 

IF   (C)  31,32,32 

0026 

31 

ESTX=K*PI-ESTX 

0027 

A=K*PI/2 

0028 

B=A4PI 

0029 

GO  TO  40 

0030 

32 

ESTX=( K-l »*PI-ESTX 

0031 

B=(K-1 J*PI/2 

0032 

IF   (K-2»  33,33,34 

0033 

33 

A  =  0 

0034 

GO  TO  40 

0035 

34 

A=B-PI 

0036 

40 

X1=ESTX 

0037 

50 

DELTA=X1* (TAN( XII -C*X1 ) / (XI /COS ( XI ) **2- 

0038 

X2=X1-DELTA 

0039 

NST0P=>NST0P*1 

0040 

IF    (X2.LT.A)   GO  TO  60 

0041 

IF   (X2.GE.BJ   GO  TO  60 

0042 

IF    (ABS ( OELTA) .LE. 0.0001 )  GO  TO  80 

0043 

X1=X2 

0044 

GO  TO  70 

0045 

60 

IF   (C.GT.O)  ESTX=(ESTX+B)/2 

0046 

IF    (C.LE.O)  ESTX=(ESTX*A)/2 

0047 

X1=ESTX 

0048 

70 

IF    (NST0P.EQ.75I   GO  TO  90 

0049 

GO  TO  50 

0050 

80 

FINDX=X2 

0051 

RETURN 

0052 

90 

FINDX=XO 

0053 

WRITE(6, 102) 

0054 

102 

FOR  NAT ( 1H   ,»FINDX  =  XO«/) 

0055 

RETURN 

0056 

END 

30 
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